Nucleation-induced transition to collective motion in active systems 
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While the existence of polar ordered states in active systems is well established, the dynamics of 
the self-assembly processes are still elusive. We study a lattice gas model of self-propelled elongated 
particles interacting through excluded volume and alignment interactions, which shows a phase 
transition from an isotropic to a polar ordered state. By analyzing the ordering process we find that 
the transition is driven by the formation of a critical nucleation cluster and a subsequent coarsening 
process. Moreover, the time to establish a polar ordered state shows a power-law divergence. 
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Understanding collective motion in driven or self- 
propelled particle systems is a topic of recent interdis- 
ciplinary interest [IH3]. Coherently moving groups have 
been observed over a broad range of scales, spanning from 
micrometer-sized systems [4HTU] over millimeter-large 
granules [TT1IT3] to large cooperative animal groups [T4] . 
The ubiquity of this phenomenon raises the questions 
of how these coherently moving and ordered clusters 
arise. Coarse-grained particle-based models [15-18] 1231 
[26] have shown that interactions favoring the alignment 
of the particles' direction of motion is sufficient for the 
emergence of large scale order. As possible origins for 
the alignment, excluded volume interactions [19 , dis- 
sipative collisions [20] and hydrodynamic interactions 
(2TJ [22] have been discussed. Complementary to these 
coarse-grained agent-based models, hydrodynamic mod- 
els have been derived from mesoscopic collision rules 
[27H3T] , or by means of microscopic interactions [32H36] . 
These approaches allow to determine the kinetic coeffi- 
cients, and thereby to analyze pattern forming instabil- 
ities. Furthermore, hydrodynamic equations based on 
symmetry arguments [37H42] were derived, enabling to 
show that collections of self-propelled particles exhibit 
a true, long-range ordered, spontaneously broken sym- 
metry state. In addition, scaling exponents have been 
calculated analytically [37] [38] , and validated by agent- 
based simulations [39] . Moreover, these studies revealed 
the existence of giant number fluctuations and long-lived 
density correlations, which have also been found experi- 
mentally [10] [nj[43]. 

While previous approaches on propelled particle sys- 
tems focussed on the long-time dynamics, the assem- 
bly processes leading to collective motion remain elusive. 
To close this gap, we present an agent-based model of 
(self-)propelled particles interacting via an effective ex- 
cluded volume interaction and a local polar alignment 
field. In the absence of excluded volume interactions, 
we find a phase transition to collective motion above a 
critical packing fraction. Strong excluded volume inter- 
actions increase this critical packing fraction or even in- 
hibit the development of collective motion. Analyzing 



the dynamics of the pattern forming processes we find 
that collective motion is accompanied by a gain in free 
volume. Close to the critical packing fraction, the onset 
of collective motion requires the spontaneous formation 
of a cluster of sufficiently large mass that acts as a nucleus 
and triggers the transition to collective motion. The cor- 
responding lag-time shows a power-law divergence upon 
approaching the critical packing fraction. 

To follow the self-assembly processes over long time 
periods and for a large number of particles we consider a 
lattice gas automaton; cf. Fig. [I] Each of the N parti- 
cles is taken as a filament of length L moving at constant 
speed and with a direction performing a persistent 
random walk. To emulate a spatially isotropic system as 
close as possible we choose a hexagonal lattice, and use 
a simulation box of hexagonal shape to avoid possible 
artifacts due to symmetry breaking with the local hexag- 
onal lattice structure; all length are measured in units of 
the lattice constant. The time scale is fixed by allowing 
each filament's head to move one lattice unit per time 
step. The particles' dynamics are assumed to be fully 
determined by its head, while the tail strictly follows the 
head's trail. The persistent random walk is implemented 
by a stochastic process, where each filament's head xo 
moves to one of its five neighboring sites x a according to 
a fixed set of conditional probabilities P prw (x a |xo). For 
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FIG. 1. (color online) Illustration of the lattice gas model. 
Filaments of length L move on a hexagonal lattice with con- 
stant speed performing a persistent random walk. 
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a fixed filament length L = 10 we choose the conditional 
probabilities such that each filament performs a persis- 
tent random walk with constant unit speed and a kinetic 
persistence length £® = 8.7 comparable to L; for more 
details see the Supplementary Material. 

Our model includes an effective excluded volume in- 
teraction and a local alignment field. Both enter into 
the model by supplementing the probabilities for the un- 
perturbed persistent random walk by appropriate fac- 
tors. In contrast to recent off-lattice simulations of over- 
damped rods, solely interacting by means of excluded vol- 
ume [19 , this allows us to separately tune excluded vol- 
ume and alignment interaction and thereby study their 
interplay. For a lattice site that is already /c-times oc- 
cupied, the probability for further occupation is reduced 
by a Boltzmann factor, e _fce , where e characterizes the 
penalty for multiple occupations. Formally, for e — >> oo, 
the limit of strict excluded volume is obtained [34 . The 
limit of weak penalties for multiple occupations (low e) is 
appropriate for motility assay experiment, where filament 
crossings occur frequently [5, 6, 8, 9 . Moreover, these 
experiments also indicate that there is a local alignment 
interaction between filaments [SJI6JI8 , 9 . To emulate such 
an interaction, each particle is assigned an alignment field 
u(x, t) at its occupied and neighboring lattice sites (see 
arrows in Fig.l); it is directed along the particles con- 
tour, and is of unit length. Overlapping alignment fields 
of different particles are averaged. In a collision event 
the alignment field modifies the transition probability to 
move from xq to x by the Boltzmann factor, a cos ^, where 
if is the relative angle between the alignment field and 
the direction of motion of the respective collision partner, 
ip = Z(u(x),x — xq). The parameter a > 1 characterizes 
the strength of the polar alignment, where stronger inter- 
actions correspond to increasing values of a. Assuming 
that all these contributions to the filament's dynamics 
are statistically independent, we arrive at the following 
update rule: Given a configuration of N filaments we 
use random sequential updating and move a chosen head 
from position xq to a target position x with probability 
P(x|x ; u, k) oc P prw (x|x ) • e~ ke ■ a cos ^ u \ 

We have performed extensive numerical simulations for 
a broad range of reduced densities p = NL/A : where A is 
the area of the simulation hexagon. For low densities an 
initial statistically homogeneous distribution of filaments 
remains isotropic both in filament position and orienta- 
tion. However, for densities above some threshold, an 
initially disordered state evolves into a state with small 
coherently moving clusters that spread perpendicular to 
their direction of motion and form bands (see Fig. [2^i- 
c, and videos in Supplementary Material). The clusters 
move coherently embedded in an isotropic background of 
randomly oriented particles. These observations resem- 
ble similar results found previously in Vicsek-like mod- 
els [17] and other agent-based lattice gas models p~6j [25] . 

The pattern formation process can readily be quan- 
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FIG. 2. (color online) (a) The boundary for our lattice gas 
model is a reflecting hexagon with a side length Lbox = 1000. 
The circle indicates the ROI with a radius of Lbox/10. Start- 
ing from a disordered state first small clusters emerge (b) 
which subsequently evolve into bands (c) (see also videos in 
Supplementary Material). Directions of cluster movement are 
indicated by arrows, (d) The polarization T^roi [orange, light 
grey] is calculated over all particles in the ROI, while the po- 
larization "Psim [blue, dark grey] is determined over all parti- 
cles in the entire simulation box. The bold and dashed arrows 
correspond to snapshots (b) and (c), respectively, (e) Orien- 
tational correlation function C(r;t) as a function of time r for 
a series of initial times t as indicated in the graph. Simulation 
results refer to the parameters a = 5, e = 0.5, and p — 0.6. 
All simulations were typically run for ~ 2 • 10 4 time steps. 



tilled by introducing the time-dependent mean polarity 
V(t), i.e. the magnitude of the polarization averaged 
over all single particle polarizations n^(t) contained in 
a region of interest (ROI) or the entire simulation box 
(SIM), respectively. The polarization of the i-th fila- 
ment is defined as the difference vector of the position 
of the filament's head between two successive time steps: 
Ui(t) := Xi(t) — Xi(t — 1) . These two measures of polar- 
ity provide complementary information about the pat- 
tern forming process, cf. Fig. [2]i. Vroi illustrates the 
formation of small polarized clusters that develop into 
bands: the sequence of spikes in the time traces indicate 
polar-moving clusters passing through the ROI In con- 
trast, Vsim characterizes the evolution of global order in 
the system: the mean value slowly grows in time and 
fluctuations are much less pronounced. These latter fluc- 
tuations originate from changes in relative polarization 
of just a few large polar-moving structures with time. 

The temporal evolution of coherent motion can be ac- 
cessed by quantifying the directional correlations in the 
filament's co- moving frame by means of the Lagrangian 
orientational correlation function: 



C(r;t) := 



(1) 
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FIG. 3. (color online) (a) Time traces of r for p = 0.6, and a set of values: e = 0, 0.5, 1 with a = 5, and e = 1 with a = 1, from 
bottom to top. Close to the phase boundary from the isotropic (ISO) to the polar ordered state (PO) [e = 1; (red)], a lag-phase 
exists, where the systems waits for the nucleation of a cluster that triggers the emergence of polar order; T nuc and T sta t are 
indicated by vertical dotted or dashed lines. Inset: for a = 1, no order develops, irrespective of the value of e (e = 0,5, 100 
from bottom to top), (b) Free volume gain $ = p — T (triangles) for a — 5 and e = 0.5 exhibits a "jump" at p c ~ 0.4 (dashed 
line) and then is significantly above the reference curve p — Thorn (p) (black + are simulation results). T sta t diverges for densities 
slightly above p c (squares), (c) Phase diagram as a function of p and e for a = 5 (solid) and a = 10 (dashed). 



The sum extends over all particles in the system, and 
t denotes the waiting-time passed since the preparation 
of the system in a state with a uniform distribution in 
both, position and orientation. Starting from a disor- 
dered state at t = 0, the directional correlation function 
decays exponentially in time r with a decay length cor- 
responding to the persistence length of the undisturbed 
persistent random walk f®. As we increase the waiting- 
time t a clear shoulder in C(r; t) builds up until it finally 
becomes stationary; for the parameters used in Fig. [2^ 
this occurs at t « 5000. Then, in the stationary regime, 
there are clearly two decay processes with well separated 
time scales: an initial fast decay and an extended time 
window with a significantly slower but again exponential 
decay. This indicates a phase separation where a fraction 
of the filaments is in a low density isotropic phase and the 
remainder is organized in coherently moving polar clus- 
ters. Filaments in the low density phase are responsible 
for the fast initial decay which has the same slope as the 
initial decay starting from an initially fully disordered 
filament configuration. In contrast, orientational corre- 
lations of filaments in clusters are significantly enhanced 
since they preferentially move into the same direction 
as their neighbors due to the impact of the alignment 
field. This results in an increased decorrelation time of 
C. Since filaments are moving at unit velocity this time 
scale can also be interpreted as a persistence length of 
the trajectories traced out by the filaments, ££ o11 . We 
attribute the loss of correlations for filaments in clusters 
mainly to reorientations of the clusters as a whole and 
to a lesser extent to collisions between clusters. The lat- 
ter and collisions of clusters with the confining boundary 
of the simulation box are responsible for the steep drop 
in the correlation function at very large times. The col- 
lective persistence length of the clusters ££ o11 depend on 
the density p as well as the strength of excluded volume 



interaction e; see Supplementary Material. 

To further analyze the self-assembly process and to 
characterize the ensuing stationary state we consider the 
fraction of occupied lattice sites T, or equivalently the 
free volume fraction := p — T. In the absence of polar 
alignment processes between the filaments, we find that 
T fluctuates around a constant value depending on the 
strength of the excluded volume interaction e (see Fig. [3^i 
inset). For a non-interacting system (e = 0, a = 1), the 
stationary value of T as a function of reduced density p 
is given by I^om (p) = 1 — e _p ; see Fig. [3]o. This result 
has previously been obtained in continuum percolation 
theory for the fraction of lattice sites occupied by statis- 
tically homogeneous distribution of overlapping rods [44] . 
In the presence of excluded volume interaction, the actual 
value of T is only slightly smaller than Thom and shows 
the same functional form. Adding polar alignment inter- 
actions leads to qualitatively different behavior. As can 
be inferred from Fig. [3)3, there is now a threshold density 
p c , where the free volume gain jumps to a value much 
larger than the corresponding value for a non-interacting 
system. Moreover, measuring the polar order parame- 
ters it turns out that this jump actually coincides with 
the onset of polar order, cf. Fig. [2^. The jump in may, 
therefore, be taken as a signature to map out the phase 
diagram as shown in Fig. [3^; for a more detailed dis- 
cussion see the Supplementary Material. It also clearly 
indicates that filaments must have formed some clusters 
much denser than expected for a purely statistical over- 
lap of filaments. Figure [3^i shows that clusters form by 
a nucleation process: There is a lag-phase during which 
T stays high close to a value obtained in the absence of 
polar alignment processes. Subsequently, there is a time 
window [T nuc , T stat ], where the available free volume frac- 
tion p — T increases towards a higher stationary value, cf. 
Fig. [3J3- The time to reach the final polar steady state 
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FIG. 4. (color online) (a) Cluster mass distribution P(m) for p — 0.6 during the lag-phase and in the stationary regime of the 
polar ordered state sampled over more than 10 4 realizations. For reference, power law distributions P(m) oc m~ s with £ = 4,5 
are indicated by dashed lines, (b) Time traces of m max [red solid], (m) [green dashed] and Sm [blue dotted] of the P(m) for 
p = 0.6. (c) Mean values of the largest cluster, (M max ) [red squares] and (M nuc ) [blue dots], serve as lower and upper bounds 
for critical nucleation cluster mass. In the transition window (shaded grey area) from the ISO to the PO phase, the fraction of 
realizations S which develop polar oder (solid line) increases from to 1. Parameters for all graphs are a = 5 and e = 0.5. 



diverges as a power law T stat ex (p — p c ) with 1. 

To further scrutinize the nucleation processes, we an- 
alyzed the cluster mass distribution P(m) and the size 
of the critical nucleation cluster as a function of den- 
sity, cf. Fig. [4]3 and c. Here clusters are defined as con- 
nected areas on the lattice with local occupation numbers 
k > k c = 5, see Supplementary Material. For an illustra- 
tion, Figure (4^i shows time traces for the maximum of the 
cluster mass distribution m max , the mean (m) and stan- 
dard deviation Sm = y/ (m 2 ) — (m) 2 close to the phase 
transition from an isotropic to an ordered state. We ob- 
serve that during the lag-phase clusters are continually 
created and destroyed such that all cluster characteris- 
tics remain constant on average and there is no gain in 
available free volume. The corresponding distribution of 
cluster sizes P(m) shows power law tails with an expo- 
nent close to S = 5, cf. Figure (4|d. For the representative 
trace shown in Fig. [4^i, there is a spontaneous nucleation 
of a large cluster at t « 2800, that albeit it partially dis- 
assembles still seems to trigger the transition to cluster 
formation. It defines the onset of a steady increase in 
(m) and also a change in the power law of the cluster 
size distribution, cf. Figure ^p. In conclusion, close to 
the critical point the system is waiting for a large enough 
density fluctuation that acts as a nucleus for cluster for- 
mation, which in turn triggers the transition to the polar 
ordered state; as shown above the waiting time diverges 
upon approaching the transition, cf. Fig. [3J3- 

To quantify the size of the critical nucleation cluster 
and further substantiate that the nucleation phenomenon 
is the mechanism driving the transition we have per- 
formed a detailed analysis of the largest cluster formed 
during the nucleation process and thereby quantified the 
heuristic observations made in Fig. [4^i. In the isotropic 
phase, one expects that the largest clusters formed over 
the whole observation window M max := max ra max (t) 

are typically too small to act as a nucleus for the emer- 



gence of polar order. However, in a finite system, there 
will still be rare occasions where fluctuations are strong 
enough to form a critical nucleus. The fraction S of re- 
alizations which develop polar order sharply increases as 
the density is raised throughout the transition window 
(see Fig. [4^), with the critical density p c « 0.4. Below 
the critical density, the values of M max for those realiza- 
tions which did not lead to polar order are a measure 
of cluster sizes which are too small to trigger the pat- 
tern forming process, and can thereby be used as lower 
bounds for the size of the critical nucleus. Figure [4^ 
shows the average value of M max for those realizations 
which did not lead to polar order. We observe that this 
quantity increases upon approaching the critical density 
from below. Above the critical density, most realizations 
will develop polar order and one may estimate the criti- 
cal nucleus size as the mean size of the largest cluster at 
the nucleation time M nuc := m max (t = T nuc ). The corre- 
sponding average value, (M nuc ) gives the typical cluster 
mass that has been large enough to trigger a transition 
to polar order, and hence serves as an upper bound for 
the critical nucleus. 

While previous numerical studies on propelled particle 
systems focussed on the long-time dynamics, we investi- 
gated the dynamics of the self-assembly process leading 
to the emergence of polar ordered states. We find that 
the ordering process is driven by the formation of a crit- 
ical nucleation cluster: once a cluster of sufficiently large 
mass has assembled, the entire system builds up collec- 
tive motion. This is reflected in the existence of a lag 
time, where the system does not succeed to assemble a 
cluster of sufficient mass. The signatures of the ensuing 
phase transitions are a jump in the free volume and a 
power-law divergence in the characteristic time for the 
formation of the polar ordered state. We suggest that 
exploring the principles underlying the self-assembly pro- 
cess of polar clusters is a promising route to understand 
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ordering phenomena in a broader class of active systems. 
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